#include "bff/project/Distortion.h"
#include <limits>

namespace bff {

std::vector<FaceData<double>> Distortion::distortion;

double computeQuasiConformalError(std::vector<Vector>& p, std::vector<Vector>& q)
{
	// compute edge vectors
	Vector u1 = p[1] - p[0];
	Vector u2 = p[2] - p[0];

	Vector v1 = q[1] - q[0];
	Vector v2 = q[2] - q[0];

	// compute orthonormal bases
	Vector e1 = u1; e1.normalize();
	Vector e2 = (u2 - dot(u2, e1)*e1); e2.normalize();

	Vector f1 = v1; f1.normalize();
	Vector f2 = (v2 - dot(v2, f1)*f1); f2.normalize();

	// project onto bases
	p[0] = Vector(0, 0, 0);
	p[1] = Vector(dot(u1, e1), dot(u1, e2), 0);
	p[2] = Vector(dot(u2, e1), dot(u2, e2), 0);

	q[0] = Vector(0, 0, 0);
	q[1] = Vector(dot(v1, f1), dot(v1, f2), 0);
	q[2] = Vector(dot(v2, f1), dot(v2, f2), 0);

	double A = 2.0*cross(u1, u2).norm();

	Vector Ss = (q[0]*(p[1].y - p[2].y) + q[1]*(p[2].y - p[0].y) + q[2]*(p[0].y - p[1].y)) / A;
	Vector St = (q[0]*(p[2].x - p[1].x) + q[1]*(p[0].x - p[2].x) + q[2]*(p[1].x - p[0].x)) / A;
	double a = dot(Ss, Ss);
	double b = dot(Ss, St);
	double c = dot(St, St);
	double det = std::sqrt(std::pow(a - c, 2) + 4.0*b*b);
	double Gamma = std::sqrt(0.5*(a + c + det));
	double gamma = std::sqrt(0.5*(a + c - det));

	if (Gamma < gamma) std::swap(Gamma, gamma);

	return Gamma/gamma;
}

Vector hsv(double h, double s, double v)
{
	double r = 0, g = 0, b = 0;

	if (s == 0) {
		r = v;
		g = v;
		b = v;

	} else {
		h = (h == 1 ? 0 : h) * 6;

		int i = (int)std::floor(h);

		double f = h - i;
		double p = v*(1 - s);
		double q = v*(1 - (s*f));
		double t = v*(1 - s*(1 - f));

		switch (i) {
			case 0:
				r = v;
				g = t;
				b = p;
				break;

			case 1:
				r = q;
				g = v;
				b = p;
				break;

			case 2:
				r = p;
				g = v;
				b = t;
				break;

			case 3:
				r = p;
				g = q;
				b = v;
				break;

			case 4:
				r = t;
				g = p;
				b = v;
				break;

			case 5:
				r = v;
				g = p;
				b = q;
				break;

			default:
				break;
		}
	}

	return Vector(r, g, b);
}

Vector Distortion::computeQuasiConformalError(const Model& model)
{
	double totalQc = 0.0;
	double minQc = std::numeric_limits<double>::infinity();
	double maxQc = -std::numeric_limits<double>::infinity();
	double totalArea = 1e-8;

	distortion.clear();
	distortion.reserve(model.size());

	for (int i = 0; i < model.size(); i++) {
		const std::vector<Face>& faces = model[i].faces;
		distortion.emplace_back(FaceData<double>(model[i]));
		std::vector<Vector> p(3), q(3);

		for (FaceCIter f = faces.begin(); f != faces.end(); f++) {
			int j = 0;
			HalfEdgeCIter he = f->halfEdge();
			do {
				p[j] = he->vertex()->position;
				q[j] = he->next()->corner()->uv;
				j++;

				he = he->next();
			} while(he != f->halfEdge());

			double qc = bff::computeQuasiConformalError(p, q);
			double a = area(f);
			totalQc += qc*a;
			totalArea += a;

			// clamp distortion to range [1, 1.5]
			distortion[i][f] = std::max(1.0, std::min(1.5, qc));
			maxQc = std::max(maxQc, qc);
			minQc = std::min(minQc, qc);
		}
	}

	return Vector(minQc, maxQc, totalQc/totalArea);
}

double computeAreaScaling(const std::vector<Vector>& p, const std::vector<Vector>& q)
{
	// Compute edge vectors and areas
	Vector u1 = p[1] - p[0];
	Vector u2 = p[2] - p[0];
	double Area = cross(u1, u2).norm();

	Vector v1 = q[1] - q[0];
	Vector v2 = q[2] - q[0];
	double area = cross(v1, v2).norm();

	return std::log(area/Area);
}

Vector seismic(double u)
{
	static std::vector<Vector> values = {{0.000, 0.000, 0.300},
									{0.000, 0.000, 0.300},
									{0.000, 0.000, 0.311},
									{0.000, 0.000, 0.311},
									{0.000, 0.000, 0.322},
									{0.000, 0.000, 0.322},
									{0.000, 0.000, 0.333},
									{0.000, 0.000, 0.333},
									{0.000, 0.000, 0.344},
									{0.000, 0.000, 0.344},
									{0.000, 0.000, 0.355},
									{0.000, 0.000, 0.355},
									{0.000, 0.000, 0.366},
									{0.000, 0.000, 0.366},
									{0.000, 0.000, 0.377},
									{0.000, 0.000, 0.377},
									{0.000, 0.000, 0.388},
									{0.000, 0.000, 0.388},
									{0.000, 0.000, 0.399},
									{0.000, 0.000, 0.399},
									{0.000, 0.000, 0.410},
									{0.000, 0.000, 0.410},
									{0.000, 0.000, 0.421},
									{0.000, 0.000, 0.421},
									{0.000, 0.000, 0.432},
									{0.000, 0.000, 0.432},
									{0.000, 0.000, 0.443},
									{0.000, 0.000, 0.443},
									{0.000, 0.000, 0.454},
									{0.000, 0.000, 0.454},
									{0.000, 0.000, 0.465},
									{0.000, 0.000, 0.465},
									{0.000, 0.000, 0.476},
									{0.000, 0.000, 0.476},
									{0.000, 0.000, 0.487},
									{0.000, 0.000, 0.487},
									{0.000, 0.000, 0.498},
									{0.000, 0.000, 0.498},
									{0.000, 0.000, 0.509},
									{0.000, 0.000, 0.509},
									{0.000, 0.000, 0.520},
									{0.000, 0.000, 0.520},
									{0.000, 0.000, 0.531},
									{0.000, 0.000, 0.531},
									{0.000, 0.000, 0.542},
									{0.000, 0.000, 0.542},
									{0.000, 0.000, 0.553},
									{0.000, 0.000, 0.553},
									{0.000, 0.000, 0.564},
									{0.000, 0.000, 0.564},
									{0.000, 0.000, 0.575},
									{0.000, 0.000, 0.575},
									{0.000, 0.000, 0.585},
									{0.000, 0.000, 0.585},
									{0.000, 0.000, 0.596},
									{0.000, 0.000, 0.596},
									{0.000, 0.000, 0.607},
									{0.000, 0.000, 0.607},
									{0.000, 0.000, 0.618},
									{0.000, 0.000, 0.618},
									{0.000, 0.000, 0.629},
									{0.000, 0.000, 0.629},
									{0.000, 0.000, 0.640},
									{0.000, 0.000, 0.640},
									{0.000, 0.000, 0.651},
									{0.000, 0.000, 0.651},
									{0.000, 0.000, 0.662},
									{0.000, 0.000, 0.662},
									{0.000, 0.000, 0.673},
									{0.000, 0.000, 0.673},
									{0.000, 0.000, 0.684},
									{0.000, 0.000, 0.684},
									{0.000, 0.000, 0.695},
									{0.000, 0.000, 0.695},
									{0.000, 0.000, 0.706},
									{0.000, 0.000, 0.706},
									{0.000, 0.000, 0.717},
									{0.000, 0.000, 0.717},
									{0.000, 0.000, 0.728},
									{0.000, 0.000, 0.728},
									{0.000, 0.000, 0.739},
									{0.000, 0.000, 0.739},
									{0.000, 0.000, 0.750},
									{0.000, 0.000, 0.750},
									{0.000, 0.000, 0.761},
									{0.000, 0.000, 0.761},
									{0.000, 0.000, 0.772},
									{0.000, 0.000, 0.772},
									{0.000, 0.000, 0.783},
									{0.000, 0.000, 0.783},
									{0.000, 0.000, 0.794},
									{0.000, 0.000, 0.794},
									{0.000, 0.000, 0.805},
									{0.000, 0.000, 0.805},
									{0.000, 0.000, 0.816},
									{0.000, 0.000, 0.816},
									{0.000, 0.000, 0.827},
									{0.000, 0.000, 0.827},
									{0.000, 0.000, 0.838},
									{0.000, 0.000, 0.838},
									{0.000, 0.000, 0.849},
									{0.000, 0.000, 0.849},
									{0.000, 0.000, 0.860},
									{0.000, 0.000, 0.860},
									{0.000, 0.000, 0.871},
									{0.000, 0.000, 0.871},
									{0.000, 0.000, 0.882},
									{0.000, 0.000, 0.882},
									{0.000, 0.000, 0.893},
									{0.000, 0.000, 0.893},
									{0.000, 0.000, 0.904},
									{0.000, 0.000, 0.904},
									{0.000, 0.000, 0.915},
									{0.000, 0.000, 0.915},
									{0.000, 0.000, 0.926},
									{0.000, 0.000, 0.926},
									{0.000, 0.000, 0.937},
									{0.000, 0.000, 0.937},
									{0.000, 0.000, 0.948},
									{0.000, 0.000, 0.948},
									{0.000, 0.000, 0.959},
									{0.000, 0.000, 0.959},
									{0.000, 0.000, 0.970},
									{0.000, 0.000, 0.970},
									{0.000, 0.000, 0.981},
									{0.000, 0.000, 0.981},
									{0.000, 0.000, 0.992},
									{0.000, 0.000, 0.992},
									{0.004, 0.004, 1.000},
									{0.004, 0.004, 1.000},
									{0.020, 0.020, 1.000},
									{0.020, 0.020, 1.000},
									{0.035, 0.035, 1.000},
									{0.035, 0.035, 1.000},
									{0.051, 0.051, 1.000},
									{0.051, 0.051, 1.000},
									{0.067, 0.067, 1.000},
									{0.067, 0.067, 1.000},
									{0.082, 0.082, 1.000},
									{0.082, 0.082, 1.000},
									{0.098, 0.098, 1.000},
									{0.098, 0.098, 1.000},
									{0.114, 0.114, 1.000},
									{0.114, 0.114, 1.000},
									{0.129, 0.129, 1.000},
									{0.129, 0.129, 1.000},
									{0.145, 0.145, 1.000},
									{0.145, 0.145, 1.000},
									{0.161, 0.161, 1.000},
									{0.161, 0.161, 1.000},
									{0.176, 0.176, 1.000},
									{0.176, 0.176, 1.000},
									{0.192, 0.192, 1.000},
									{0.192, 0.192, 1.000},
									{0.208, 0.208, 1.000},
									{0.208, 0.208, 1.000},
									{0.224, 0.224, 1.000},
									{0.224, 0.224, 1.000},
									{0.239, 0.239, 1.000},
									{0.239, 0.239, 1.000},
									{0.255, 0.255, 1.000},
									{0.255, 0.255, 1.000},
									{0.271, 0.271, 1.000},
									{0.271, 0.271, 1.000},
									{0.286, 0.286, 1.000},
									{0.286, 0.286, 1.000},
									{0.302, 0.302, 1.000},
									{0.302, 0.302, 1.000},
									{0.318, 0.318, 1.000},
									{0.318, 0.318, 1.000},
									{0.333, 0.333, 1.000},
									{0.333, 0.333, 1.000},
									{0.349, 0.349, 1.000},
									{0.349, 0.349, 1.000},
									{0.365, 0.365, 1.000},
									{0.365, 0.365, 1.000},
									{0.380, 0.380, 1.000},
									{0.380, 0.380, 1.000},
									{0.396, 0.396, 1.000},
									{0.396, 0.396, 1.000},
									{0.412, 0.412, 1.000},
									{0.412, 0.412, 1.000},
									{0.427, 0.427, 1.000},
									{0.427, 0.427, 1.000},
									{0.443, 0.443, 1.000},
									{0.443, 0.443, 1.000},
									{0.459, 0.459, 1.000},
									{0.459, 0.459, 1.000},
									{0.475, 0.475, 1.000},
									{0.475, 0.475, 1.000},
									{0.490, 0.490, 1.000},
									{0.490, 0.490, 1.000},
									{0.506, 0.506, 1.000},
									{0.506, 0.506, 1.000},
									{0.522, 0.522, 1.000},
									{0.522, 0.522, 1.000},
									{0.537, 0.537, 1.000},
									{0.537, 0.537, 1.000},
									{0.553, 0.553, 1.000},
									{0.553, 0.553, 1.000},
									{0.569, 0.569, 1.000},
									{0.569, 0.569, 1.000},
									{0.584, 0.584, 1.000},
									{0.584, 0.584, 1.000},
									{0.600, 0.600, 1.000},
									{0.600, 0.600, 1.000},
									{0.616, 0.616, 1.000},
									{0.616, 0.616, 1.000},
									{0.631, 0.631, 1.000},
									{0.631, 0.631, 1.000},
									{0.647, 0.647, 1.000},
									{0.647, 0.647, 1.000},
									{0.663, 0.663, 1.000},
									{0.663, 0.663, 1.000},
									{0.678, 0.678, 1.000},
									{0.678, 0.678, 1.000},
									{0.694, 0.694, 1.000},
									{0.694, 0.694, 1.000},
									{0.710, 0.710, 1.000},
									{0.710, 0.710, 1.000},
									{0.725, 0.725, 1.000},
									{0.725, 0.725, 1.000},
									{0.741, 0.741, 1.000},
									{0.741, 0.741, 1.000},
									{0.757, 0.757, 1.000},
									{0.757, 0.757, 1.000},
									{0.773, 0.773, 1.000},
									{0.773, 0.773, 1.000},
									{0.788, 0.788, 1.000},
									{0.788, 0.788, 1.000},
									{0.804, 0.804, 1.000},
									{0.804, 0.804, 1.000},
									{0.820, 0.820, 1.000},
									{0.820, 0.820, 1.000},
									{0.835, 0.835, 1.000},
									{0.835, 0.835, 1.000},
									{0.851, 0.851, 1.000},
									{0.851, 0.851, 1.000},
									{0.867, 0.867, 1.000},
									{0.867, 0.867, 1.000},
									{0.882, 0.882, 1.000},
									{0.882, 0.882, 1.000},
									{0.898, 0.898, 1.000},
									{0.898, 0.898, 1.000},
									{0.914, 0.914, 1.000},
									{0.914, 0.914, 1.000},
									{0.929, 0.929, 1.000},
									{0.929, 0.929, 1.000},
									{0.945, 0.945, 1.000},
									{0.945, 0.945, 1.000},
									{0.961, 0.961, 1.000},
									{0.961, 0.961, 1.000},
									{0.976, 0.976, 1.000},
									{0.976, 0.976, 1.000},
									{0.992, 0.992, 1.000},
									{0.992, 0.992, 1.000},
									{1.000, 0.992, 0.992},
									{1.000, 0.992, 0.992},
									{1.000, 0.976, 0.976},
									{1.000, 0.976, 0.976},
									{1.000, 0.961, 0.961},
									{1.000, 0.961, 0.961},
									{1.000, 0.945, 0.945},
									{1.000, 0.945, 0.945},
									{1.000, 0.929, 0.929},
									{1.000, 0.929, 0.929},
									{1.000, 0.914, 0.914},
									{1.000, 0.914, 0.914},
									{1.000, 0.898, 0.898},
									{1.000, 0.898, 0.898},
									{1.000, 0.882, 0.882},
									{1.000, 0.882, 0.882},
									{1.000, 0.867, 0.867},
									{1.000, 0.867, 0.867},
									{1.000, 0.851, 0.851},
									{1.000, 0.851, 0.851},
									{1.000, 0.835, 0.835},
									{1.000, 0.835, 0.835},
									{1.000, 0.820, 0.820},
									{1.000, 0.820, 0.820},
									{1.000, 0.804, 0.804},
									{1.000, 0.804, 0.804},
									{1.000, 0.788, 0.788},
									{1.000, 0.788, 0.788},
									{1.000, 0.773, 0.773},
									{1.000, 0.773, 0.773},
									{1.000, 0.757, 0.757},
									{1.000, 0.757, 0.757},
									{1.000, 0.741, 0.741},
									{1.000, 0.741, 0.741},
									{1.000, 0.725, 0.725},
									{1.000, 0.725, 0.725},
									{1.000, 0.710, 0.710},
									{1.000, 0.710, 0.710},
									{1.000, 0.694, 0.694},
									{1.000, 0.694, 0.694},
									{1.000, 0.678, 0.678},
									{1.000, 0.678, 0.678},
									{1.000, 0.663, 0.663},
									{1.000, 0.663, 0.663},
									{1.000, 0.647, 0.647},
									{1.000, 0.647, 0.647},
									{1.000, 0.631, 0.631},
									{1.000, 0.631, 0.631},
									{1.000, 0.616, 0.616},
									{1.000, 0.616, 0.616},
									{1.000, 0.600, 0.600},
									{1.000, 0.600, 0.600},
									{1.000, 0.584, 0.584},
									{1.000, 0.584, 0.584},
									{1.000, 0.569, 0.569},
									{1.000, 0.569, 0.569},
									{1.000, 0.553, 0.553},
									{1.000, 0.553, 0.553},
									{1.000, 0.537, 0.537},
									{1.000, 0.537, 0.537},
									{1.000, 0.522, 0.522},
									{1.000, 0.522, 0.522},
									{1.000, 0.506, 0.506},
									{1.000, 0.506, 0.506},
									{1.000, 0.490, 0.490},
									{1.000, 0.490, 0.490},
									{1.000, 0.475, 0.475},
									{1.000, 0.475, 0.475},
									{1.000, 0.459, 0.459},
									{1.000, 0.459, 0.459},
									{1.000, 0.443, 0.443},
									{1.000, 0.443, 0.443},
									{1.000, 0.427, 0.427},
									{1.000, 0.427, 0.427},
									{1.000, 0.412, 0.412},
									{1.000, 0.412, 0.412},
									{1.000, 0.396, 0.396},
									{1.000, 0.396, 0.396},
									{1.000, 0.380, 0.380},
									{1.000, 0.380, 0.380},
									{1.000, 0.365, 0.365},
									{1.000, 0.365, 0.365},
									{1.000, 0.349, 0.349},
									{1.000, 0.349, 0.349},
									{1.000, 0.333, 0.333},
									{1.000, 0.333, 0.333},
									{1.000, 0.318, 0.318},
									{1.000, 0.318, 0.318},
									{1.000, 0.302, 0.302},
									{1.000, 0.302, 0.302},
									{1.000, 0.286, 0.286},
									{1.000, 0.286, 0.286},
									{1.000, 0.271, 0.271},
									{1.000, 0.271, 0.271},
									{1.000, 0.255, 0.255},
									{1.000, 0.255, 0.255},
									{1.000, 0.239, 0.239},
									{1.000, 0.239, 0.239},
									{1.000, 0.224, 0.224},
									{1.000, 0.224, 0.224},
									{1.000, 0.208, 0.208},
									{1.000, 0.208, 0.208},
									{1.000, 0.192, 0.192},
									{1.000, 0.192, 0.192},
									{1.000, 0.176, 0.176},
									{1.000, 0.176, 0.176},
									{1.000, 0.161, 0.161},
									{1.000, 0.161, 0.161},
									{1.000, 0.145, 0.145},
									{1.000, 0.145, 0.145},
									{1.000, 0.129, 0.129},
									{1.000, 0.129, 0.129},
									{1.000, 0.114, 0.114},
									{1.000, 0.114, 0.114},
									{1.000, 0.098, 0.098},
									{1.000, 0.098, 0.098},
									{1.000, 0.082, 0.082},
									{1.000, 0.082, 0.082},
									{1.000, 0.067, 0.067},
									{1.000, 0.067, 0.067},
									{1.000, 0.051, 0.051},
									{1.000, 0.051, 0.051},
									{1.000, 0.035, 0.035},
									{1.000, 0.035, 0.035},
									{1.000, 0.020, 0.020},
									{1.000, 0.020, 0.020},
									{1.000, 0.004, 0.004},
									{1.000, 0.004, 0.004},
									{0.994, 0.000, 0.000},
									{0.994, 0.000, 0.000},
									{0.986, 0.000, 0.000},
									{0.986, 0.000, 0.000},
									{0.978, 0.000, 0.000},
									{0.978, 0.000, 0.000},
									{0.971, 0.000, 0.000},
									{0.971, 0.000, 0.000},
									{0.963, 0.000, 0.000},
									{0.963, 0.000, 0.000},
									{0.955, 0.000, 0.000},
									{0.955, 0.000, 0.000},
									{0.947, 0.000, 0.000},
									{0.947, 0.000, 0.000},
									{0.939, 0.000, 0.000},
									{0.939, 0.000, 0.000},
									{0.931, 0.000, 0.000},
									{0.931, 0.000, 0.000},
									{0.924, 0.000, 0.000},
									{0.924, 0.000, 0.000},
									{0.916, 0.000, 0.000},
									{0.916, 0.000, 0.000},
									{0.908, 0.000, 0.000},
									{0.908, 0.000, 0.000},
									{0.900, 0.000, 0.000},
									{0.900, 0.000, 0.000},
									{0.892, 0.000, 0.000},
									{0.892, 0.000, 0.000},
									{0.884, 0.000, 0.000},
									{0.884, 0.000, 0.000},
									{0.876, 0.000, 0.000},
									{0.876, 0.000, 0.000},
									{0.869, 0.000, 0.000},
									{0.869, 0.000, 0.000},
									{0.861, 0.000, 0.000},
									{0.861, 0.000, 0.000},
									{0.853, 0.000, 0.000},
									{0.853, 0.000, 0.000},
									{0.845, 0.000, 0.000},
									{0.845, 0.000, 0.000},
									{0.837, 0.000, 0.000},
									{0.837, 0.000, 0.000},
									{0.829, 0.000, 0.000},
									{0.829, 0.000, 0.000},
									{0.822, 0.000, 0.000},
									{0.822, 0.000, 0.000},
									{0.814, 0.000, 0.000},
									{0.814, 0.000, 0.000},
									{0.806, 0.000, 0.000},
									{0.806, 0.000, 0.000},
									{0.798, 0.000, 0.000},
									{0.798, 0.000, 0.000},
									{0.790, 0.000, 0.000},
									{0.790, 0.000, 0.000},
									{0.782, 0.000, 0.000},
									{0.782, 0.000, 0.000},
									{0.775, 0.000, 0.000},
									{0.775, 0.000, 0.000},
									{0.767, 0.000, 0.000},
									{0.767, 0.000, 0.000},
									{0.759, 0.000, 0.000},
									{0.759, 0.000, 0.000},
									{0.751, 0.000, 0.000},
									{0.751, 0.000, 0.000},
									{0.743, 0.000, 0.000},
									{0.743, 0.000, 0.000},
									{0.735, 0.000, 0.000},
									{0.735, 0.000, 0.000},
									{0.727, 0.000, 0.000},
									{0.727, 0.000, 0.000},
									{0.720, 0.000, 0.000},
									{0.720, 0.000, 0.000},
									{0.712, 0.000, 0.000},
									{0.712, 0.000, 0.000},
									{0.704, 0.000, 0.000},
									{0.704, 0.000, 0.000},
									{0.696, 0.000, 0.000},
									{0.696, 0.000, 0.000},
									{0.688, 0.000, 0.000},
									{0.688, 0.000, 0.000},
									{0.680, 0.000, 0.000},
									{0.680, 0.000, 0.000},
									{0.673, 0.000, 0.000},
									{0.673, 0.000, 0.000},
									{0.665, 0.000, 0.000},
									{0.665, 0.000, 0.000},
									{0.657, 0.000, 0.000},
									{0.657, 0.000, 0.000},
									{0.649, 0.000, 0.000},
									{0.649, 0.000, 0.000},
									{0.641, 0.000, 0.000},
									{0.641, 0.000, 0.000},
									{0.633, 0.000, 0.000},
									{0.633, 0.000, 0.000},
									{0.625, 0.000, 0.000},
									{0.625, 0.000, 0.000},
									{0.618, 0.000, 0.000},
									{0.618, 0.000, 0.000},
									{0.610, 0.000, 0.000},
									{0.610, 0.000, 0.000},
									{0.602, 0.000, 0.000},
									{0.602, 0.000, 0.000},
									{0.594, 0.000, 0.000},
									{0.594, 0.000, 0.000},
									{0.586, 0.000, 0.000},
									{0.586, 0.000, 0.000},
									{0.578, 0.000, 0.000},
									{0.578, 0.000, 0.000},
									{0.571, 0.000, 0.000},
									{0.571, 0.000, 0.000},
									{0.563, 0.000, 0.000},
									{0.563, 0.000, 0.000},
									{0.555, 0.000, 0.000},
									{0.555, 0.000, 0.000},
									{0.547, 0.000, 0.000},
									{0.547, 0.000, 0.000},
									{0.539, 0.000, 0.000},
									{0.539, 0.000, 0.000},
									{0.531, 0.000, 0.000},
									{0.531, 0.000, 0.000},
									{0.524, 0.000, 0.000},
									{0.524, 0.000, 0.000},
									{0.516, 0.000, 0.000},
									{0.516, 0.000, 0.000},
									{0.508, 0.000, 0.000},
									{0.508, 0.000, 0.000},
									{0.500, 0.000, 0.000},
									{0.500, 0.000, 0.000}};

	// find the two nearest indices in to the colormap lookup table, then
	// return a linear blend between them
	double scaled = u*(values.size() - 1);
	double lower = std::floor(scaled);
	double upperBlend = scaled - lower;
	unsigned int lowerInd = static_cast<unsigned int>(lower);
	unsigned int upperInd = lowerInd + 1;

	return (1.0 - upperBlend)*values[lowerInd] + upperBlend*values[upperInd];
}

Vector Distortion::computeAreaScaling(const std::vector<Face>& faces)
{
	double totalU = 0.0;
	double minU = std::numeric_limits<double>::infinity();
	double maxU = -std::numeric_limits<double>::infinity();
	double totalArea = 1e-8;
	std::vector<Vector> p(3), q(3);

	for (FaceCIter f = faces.begin(); f != faces.end(); f++) {
		int i = 0;
		HalfEdgeCIter he = f->halfEdge();
		do {
			p[i] = he->vertex()->position;
			q[i] = he->next()->corner()->uv;
			i++;

			he = he->next();
		} while(he != f->halfEdge());

		double u = bff::computeAreaScaling(p, q);
		maxU = std::max(maxU, u);
		minU = std::min(minU, u);
		double a = area(f);
		totalU += u*a;
		totalArea += a;
	}

	return Vector(minU, maxU, totalU/totalArea);
}

Vector Distortion::computeAreaScaling(const Model& model)
{
	double totalU = 0.0;
	double minU = std::numeric_limits<double>::infinity();
	double maxU = -std::numeric_limits<double>::infinity();
	double totalArea = 1e-8;

	distortion.clear();
	distortion.reserve(model.size());

	for (int i = 0; i < model.size(); i++) {
		const std::vector<Face>& faces = model[i].faces;
		distortion.emplace_back(FaceData<double>(model[i]));
		std::vector<Vector> p(3), q(3);

		for (FaceCIter f = faces.begin(); f != faces.end(); f++) {
			int j = 0;
			HalfEdgeCIter he = f->halfEdge();
			do {
				p[j] = he->vertex()->position;
				q[j] = he->next()->corner()->uv;
				j++;

				he = he->next();
			} while(he != f->halfEdge());

			double u = bff::computeAreaScaling(p, q);
			maxU = std::max(maxU, u);
			minU = std::min(minU, u);
			double a = area(f);
			totalU += u*a;
			totalArea += a;

			distortion[i][f] = u;
		}

		for (FaceCIter f = faces.begin(); f != faces.end(); f++) {
			distortion[i][f] = (distortion[i][f] + 15)/30;
			if (distortion[i][f] < 0) distortion[i][f] = 0;
			if (distortion[i][f] > 1) distortion[i][f] = 1;
		}
	}

	return Vector(minU, maxU, totalU/totalArea);
}

Vector Distortion::color(FaceCIter f, int meshIndex, bool conformalColors)
{
	if (conformalColors) return hsv((2.0 - 4.0*(distortion[meshIndex][f] - 1.0))/3.0, 0.7, 0.65);
	return seismic(distortion[meshIndex][f]);
}

} // namespace bff
